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Time-domain  Deconvolution 
Removes  the  Effects  of 
Near-held  Scatterers 


1.  INTRODUCTION 

This  report  studies  the  removal  of  interference  caused,  in  linear  propagation,  by 
neax-field  scatterers.  The  report’s  numerical  example  is  a  schematic  representa¬ 
tive  for  the  removal  of  electromagnetic  interference  caused  by  a  mast,  wing  tip, 
or  fin  that  is  in  the  near  field  of  a  ship  or  airborne  receiving  antenna  [1] . 

The  above  problem  will  be  shown  to  reduce  to  the  deconvolution  of  the  right- 
hand  side  of  Eq.  (1)  below,  which  could  be  easily  accomplished  using  Fourier  or 
Laplace  transforms.  This  report’s  goal,  therefore,  is  to  find  what  may  be  the  best 
alternative  time-domain  algorithm.  The  best  time-domain  algorithm  considered 
here  is  competitive  with  frequency- domain  algorithms  in  that  the  t-domain  algo¬ 
rithm  feeds  one-degree-smoothed  data  into  a  first-kind- Volterra-equation  solver 
[2]  that  is  second-order  accurate  and  stable,  and  for  which  Richardson  extrap¬ 
olation  yields  a  fourth-order-accurate  method.  Special  care  is  taken  to  find 
a  deconvolution  algorithm  that  accommodates  discontinuity-related  numerical 
noise  in  finite-difference-time-domain  (fdtd)  data. 

The  central  equation  for  the  preferred  algorithm  will  be  shown  to  be 

fRis)ds=  fKn(t-s)fi,,is)ds,  (1) 

Jo  Jo 

where  R  is  the  received  signal  that  results  from  the  corruption  of  the  incident 
signal  /inc  by  a  near-field  scatter.  The  identity  Eq.  (1)  follows  directly  from 
the  Duhamel  theorem  [3]  concerning  the  Heaviside-step  response  iCn  of  a  linear 
system.  The  discontinuous  kernel  will  be  computed  here  using  an  FDTD 
(finite  difference  time  domain)  method,  despite  the  noise  introduced  by  the  FDTD 
propagation  of  a  discontinuity.  Indeed,  this  approach  Eq.  (1)  was  stimulated  by 
an  earlier  paper  [4],  which  established  the  usefulness  of  FDTD  propagation  of 
discontinuities  in  linear  scattering.  This  report  and  [4]  are  complementary  in 
that  the  central  operator  in  [4]  is  a  convolution  and  the  central  operator  here  is 
the  deconvolution  of  Eq.  (1). 

Received  for  publication  27  May  1998. 
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A  prototype  problem  will  be  defined  in  Section  2.  Section  3  considers  an 
infinite  sequence  of  time-domain  algorithms,  each  of  which  could  solve  the  pro¬ 
totype  problem.  The  best  of  these  algorithms  is  fotmd  using  two  niimerical 
criteria.  Section  4  solves  the  prototype  problem  of  Section  2  numerically.  The 
conclusion  (Section  5)  describes  the  relation  of  the  present  work  to  an  earlier- 
published  paper  [4],  which  together  find  two  uses  for  the  FDTD  propagation  of 
discontinuous  functions.  The  potential  use  of  laboratory  data  also  is  discussed 
in  Section  5.  The  Appendix  studies  the  only  approximation  (tnmcation  of  su- 
perluminal  components)  that  is  made  in  the  best  time-domain  algorithm. 


2.  A  PROTOTYPE  PROBLEM 

A  prototype  problem  is  defined  here,  and  it  will  later  be  solved  numerically  in 
Section  4  as  an  example  of  a  more  general  procedure  developed  in  Section  3. 
This  prototype  problem  involves  the  realistic  parameters  of  an  existing  antenna 
[5],  which  is  sketched  in  Figure  1.  This  section  defines  the  problem,  says  what 
it  represents,  and  then  explains  how  this  report’s  analysis  applies  to  more  com¬ 
plicated  problems  that  involve  multiple  scattering. 
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Figure  1.  Scale  drawing  of  a  metal  disk  and  the  32  locations  (x) 
where  the  total  field  is  received  inside  the  computational  domain. 


There  is  initially,  for  t  <  0,  no  field  in  the  rectangular  domain  in  Figure  1. 
At  t  =  0  a  t-dependent  field  becomes  incident  uniformly  from  the  right,  with 
its  electric-field  component  being  always  perpendicular  to  the  page.  This  field 
represents  a  pulse  that  is  incident  from  a  distant  source,  known  also  as  a  soft 
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soiirce,  and  the  incident  pulse  is  therefore  added,  at  each  timestep,  to  the  total 
field  at  each  point  on  the  right-hand  edge  of  Figure  1.  The  incident  field  itself 
is  never  received;  instead,  it  propagates  through  free  space  (e  =  1)  and  scatters 
from  a  conducting  (a  <  oo)  object,  which  is  drawn  as  a  disk  in  Figmre  1,  and  the 
total  <-dependent  electric  field  is  then  measured  at  each  cross-marked  location 
in  the  figure.  It  is  assumed  that  the  measurements  do  not  perturb  the  field. 
Let  Rj{t)  be  the  field  that  is  measured  at  the  location.  Then  this  report’s 
prototype  problem  is  as  follows:  Given  the  total  field  Rj{t)  at  a  single  known 
location,  knowing  that  the  field  was  produced  by  a  soft-source  plane- wave  pulse 
incident  from  a  known  direction,  and  knowing  also  the  location,  shape,  and 
composition  (e  and  a)  of  the  scatterer,  one  must  then  compute  the  time  trace 
of  the  incident  pulse. 

The  protot3rpe  problem  (above)  is  a  schematic  representative  for  the  removal 
of  electromagnetic  interference  with  a  ship  or  airborne  antenna  that  is  in  the 
near  field  of  a  mast,  wing  tip,  or  fin  [1].  In  Section  4  the  incident  field  will  be 
taken  to  be  a  1-cycle  sinusoid  with  a  5.45-GHz  carrier  frequency,  whose  free- 
space  wavelength  is  approximately  the  diameter  of  the  steel  disk  in  Figure  1. 
Because  the  wavelength  is  also  comparable  to  the  length  of  two  cross-marked 
intervals  in  Figure  1,  this  prototype  problem  represents  a  32-element  phased- 
array  radar  antenna  [5]  that  has  a  wavelength-sized  metal  pipe  located  40-cm 
in  front  of  it.  That  pipe  introduces  a  reasonably  large  amount  of  interference, 
which  we  seek  to  remove. 

One  of  the  prototype  problem’s  more  idealized  assumptions  is  that  the 
antenna  measures  the  field  nonperturbatively.  This  assumption  neglects  the 
important  practical  effect  of  multiple  scattering  among  the  32  elements,  and 
it  thereby  conforms  to  the  tradition  that  prototype  problems  be  simple.  One 
way  to  simulate  multiple  scattering,  however,  would  be  to  have  an  imperfect 
absorber  in  a  small  area  near  each  x -marked  element,  and  perhaps  to  have  a 
small  conductor  behind  each  imperfect  absorber.  In  any  case,  the  fdtd  response 
would  still  be  a  linear  operator  with  respect  to  the  incident  field,  regardless  of 
what  reasonable  conductors  and  absorbers  are  used.  The  analysis  in  this  report 
wotold  a,pply  without  change  to  any  such  reasonable  linear  system,  including  all 
linear  systems  that  have  multiple  scattering. 
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3.  NUMERICAL  METHODS  AND  THEIR  PROPERTIES 


Because  the  deconvolution  problem  (Section  2)  involves  a  linear  system,  it  can 
be  easily  solved  using  Fourier  or  Laplace  transforms.  This  section  will  consider 
infinitely  many  time-domain  alternatives  (Eqs.  (2))  and  will  use  two  numerical 
criteria  to  select  the  best  alternative. 

Linearity  yields  infinitely  many  integral  equations  for  the  incident  field  /inc- 
The  equations  axe 


m 

dr'm 


=  /  Ki(t  -  i)/inc(s)<is 
Jo 

=  f  5)/inc(5)ds 

Jo 

=  [  Km{t  -  s)fijic{s)ds, 
Jo 


(2a) 

(2b) 

(2c) 


where  is  the  square  of  the  antiderivative  operator  d^^R{t)  =  J’Qi?(s)ds, 
and  where  the  subscripted  kernels  indicate  the  delta-function  response,  Ks,  the 
Heaviside  response,  Ku,  and  the  ramp-function  [tH(t)]  response.  Km,  at  the 
location  where  the  received  signal  R  is  measured.  In  linear  hyperbolic  systems, 
such  as  linear  t-domain  electromagnetics,  a  propagation-of-singularities  argu¬ 
ment  [3]  shows  that  Ks  has  a  nonzero  delta-function  component;  however,  [4] 
showed  that  a  bounded  (Loo)  approximation  to  6  could  be  propagated  usefully 
using  FDTD.  Assinning  that  Ks  in  Eq.  (2a)  would  be  computed  using  such  an 
Loo  computation,  it  follows  that  the  sequence  of  Eqs.  (2a),  (2b),  (2c),  ...  are  all 
first-kind  Volterra  integral  equations  (for  /inc)  with  convolution  kernels;  whence 
the  word  “deconvolution”  was  brought  into  the  title  of  this  report. 

We  now  consider  a  numerical  property  that  will  be  used  as  a  selection 
criterion  for  Eqs.  (2):  Linz  [2]  showed  that  if  the  left-hand  side  of  a  general 
first-kind  Volterra  equation 


L(t)  =  f  K(t,3)f{s)ds  (3) 

is  perturbed  by  an  amount  AL  then  the  resulting  perturbation  in  the  solution 
of  Eq.  (3)  is,  in  what  is  probably  the  best  case, 


A/  =  0  (h-^AL)  . 


(4) 
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That  ill-posedness  result  in  Eq.  (4)  favors  equations  that  have  smooth  left-hand 
sides;  consequently,  Eq.  (2a)  is  eliminated  as  a  candidate. 

We  still  have  infinitely  many  candidates  —  Eqs.  (2b),  (2c),  ...  —  from 
which  we  will  be  able  to  extract  the  best  candidate  only  after  a  method  for 
solving  first-kind  equations  is  described.  This  method  will  be  described  as  it 
applies  to  the  eventually  preferred  candidate  Eq.  (2b),  but  the  same  method  is 
easily  adapted  to  suit  all  candidates  in  Eqs.  (2). 

Direct  methods  for  solving  first-kind  equations  (Eq.  (3))  follow  immediately 
from  discretization  of  the  integral.  The  direct  midpoint-rule  discretization  of 
Eq.  (2b)  yields 


[A-K'h  (I)] 


n-1 


h  ^  Ku,n-ifi  5 


1  =  1 


(5a) 

(5b) 


where 


fn  =  fine  [(«  -  1)  h] 

(6a) 

^H.n  =  Ku  [(n  -  i)  h] 

(6b) 

pnh 

[a,-'fi]„=/ 

Jo 

(6c) 

are  stepsize-h  discretizations.  Notice  that  if 


\hKE{h/2)\  <C  min  (|/„|,|[ar’J2]„|) 


(7) 


then  the  parenthetical  numerator  of  Eq.  (5b)  would  be  a  small  difference  of 
large  numbers,  which  would  cause  a  loss  of  significant  digits.  Candidate  (2b) 
is  therefore  preferred  because  a  propagation-of-singularities  argument  [3]  shows 
that  iirH(0+)  ^  0,  whereas  candidates  Eqs.  (2c),  (2d),  et  seq.  have  integral 
kernels  that  axe  continuous  and  zero  at  t  =  0,  with  the  zeros  being  first  order  for 
Ktn,  second  order  for  Kfi-R,  and  so  forth.  Indeed,  numerical  experiments  have 
shown  that  the  above-described  loss  of  significant  digits  causes  a  rapid  numerical 
blowup  for  the  first-order-zero  case  (Eq.  (2c)).  Thus,  Eq.  (2b)  is  the  best  of  the 
infinitely  many  algorithms  in  Eqs.  (2). 
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We  turn  now  to  the  numerical  properties  of  the  best  method,  Eq.  (2b). 

Linz  [2]  has  shown  that  error  in  the  the  direct  midpoint-rule  solution  of 
first-kind  Volterra  equations  is 


U.«-fn  =  ^e„  +  0(h'),  (8) 

where  £„  is  independent  of  the  stepsize  h.  Examining  the  algorithm  (2b)  in  light 
of  Eqs.  (4)  and  (8),  we  see  that  smoothing  the  data  reduces  the  algorithm’s  ill- 
posedness,  and  that  the  smoothed  data  are  fed  into  an  algorithm  Eq.  (5) 

that  is  stable  and  second-order  accurate,  and  for  which  Richardson  extrapolation 
[2],  based  on  Eq.  (8),  would  yield  a  fourth-order-accurate  solution.  For  any 
given  error  criterion,  the  method  would  allow  one  to  increase  the  stepsize 
and  thereby  further  reduce  the  method’s  noise  sensitivity  according  to  Eq.  (4). 

An  obvious  challenge  in  using  Eq.  (2b)  is  the  accommodation  of  the  large 
amotmt  of  numerical  noise  that  is  caused  by  the  FDTD  propagation  of  the  discon¬ 
tinuous  function  H(t)  when  computing  Kn.  The  propagated  discontinuity  could 
have  been  avoided  by  instead  propagating  the  ramp  function  tH(t)  and  then 
differentiating,  as  in  Ku  =  dtKtR,  but  noise  would  then  have  been  introduced 
by  numerical  differentiation.  Therefore,  for  the  sake  of  simplicity,  the  kernel  ATh 
is  computed  directly  by  propagating  H(t),  which  is  similar  to  a  process  whose 
usefulness  has  already  been  established  [4] .  The  stepsize  h  is  decreased  until  the 
/i-dependent,  discontinuity-related  noise  is  contained  in  a  frequency  band  that 
is  separated  from  the  spectrum  of  the  main  physical  features  of  ATh  •  In.  practice 
(Section  4)  the  numerical  noise  consisted  of  oscillations  whose  periods  varied 
from  lOA  to  20h,  regardless  of  the  value  of  h.  A  40-point  centered  filter  $  was 
then  defined,  using  tm  =  19^  and  =  20/i  in  (Al).  The  Appendix  uses  #  to 
derive  the  identity 


(9) 

in  which  *  henceforth  represents  /I^j^-type  convolution.  Notice  that  both  Eq.  (9) 
and  its  imfiltered  predecessor  Eq.  (2b)  are  exact  for  all  t.  The  next,  and  final,  nu¬ 
merical  consideration  will,  however,  lead  to  an  approximate  version  of  Eqs.  (2b) 
and  (9). 

The  wavefront  of  an  FDTD-computed  field  travels  with  the  superliuninal 
velocity  c/cFL  because  of  the  nature  of  time  stepping  in  a  CFL-stabilized  com¬ 
putation.  The  filtered  fields  ($  *  R)  and  ($  *  A’h)  were  therefore  trimcated  as 
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{^*R) 


r 


(10a) 

(10b) 


where  t  is  now  measured  with  respect  to  the  analytically  determined  wavefront- 
ariival  time  (ta)  at  each  receiver.  (That  arrival  time  is  the  same  for  all  receivers 
in  Figure  1.  For  more  complicated  media,  the  arrival  time  can  be  computed 
using  characteristics  or  the  eikonal  equation  [3].)  The  quantities  on  the  right- 
hand  sides  of  Eqs.  (10a)  and  (10b)  are  used,  respectively,  as  the  Ke  and  R  terms 
in  Eq.  (5).  Notice  that  the  truncation  in  Eq.  (10a)  of  the  necessarily  superluminal 
FDTD-computed  Ke  assures  that  the  right-hand  side  of  Eq.  (10a)  will  have  an 
initial  discontinuity,  which  is  beneficial  from  the  standpoint  of  Eq.  (7).  This 
may  be  the  only  known  computational  advantage  of  the  superluminal  feature  of 
FDTD.  If,  however,  Ke  had  been  computed  exactly  by  analytical  means  then, 
because  of  a  propagation-of-singularities  argument  [3],  it  would  already  have 
been  discontinuous  and,  because  exact  computations  have  no  noise,  filtering 
would  not  have  been  necessary.  Regardless  of  whether  Ke  is  computed  exactly 
or  with  filtered-then- truncated  fdtd,  one  does  obtain  a  discontinuous  integral 
kernel  that,  together  with  an  exact  or  filtered- then-truncated  R,  can  be  fed  into 
Eq.  (5),  with  the  same  numerical  benefits  that  Eqs.  (4),  (7),  and  (8)  ordinarily 
yield. 

The  main  disadvantage  of  truncation  is  that  the  procedure  of  feeding  Eq.  (10) 
into  Eq.  (5)  is  inherently  approximate,  unhke  Eq.  (9)  which  uses  filtered-but- 
untruncated  data.  The  Appendix  shows  that  the  trimcation  errors  in  (10), 
however,  will  vanish  in  the  h 0  limit,  provided  that  the  underlying  FDTD 
computation  is  reasonable;  and  in  that  matter  there  is  a  fine  point  Eq.  (A16) 
that  arises  from  a  discontinuity. 

This  section  has  considered  infinitely  many  —  but  not  all  —  time-domain 
methods.  Resolvent-kernel  methods  [10],  for  example,  seem  especially  appealing 
for  processing  many  received  signals;  for  if  one  could  solve  ($  *  Ke)  *  p  =  1 
for  p,  then  (filtered)  incident  fields  could  be  quickly  inferred  using  #  *  /  = 
R*{^  *  p)  —  J?(0‘*')5t“^($  *  p).  A  numerical  experiment,  however,  showed  that 
solving  (#  *  Ke)  *  p  =  1  can  be  ill  conditioned,  but  many  other  resolvent-kernel 
methods  remain  unexplored. 


4.  NUMERICAL  SOLUTION 

This  report’s  prototype  problem  was  defined  in  Figure  1  and  Section  2.  The 
incident  field  was  a  one-cycle  5.45-GHz  sinusoidal  pulse.  That  carrier  frequency 
had  been  chosen  to  be  at  the  midpoint  of  the  operating  band  of  the  existing 
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antenna  [5]  that  this  problem  models.  The  corresponding  wavelength  was  5.5  cm, 
which  can  be  compared  in  Figure  1  with  the  3-cm  element  ( x  symbol)  spacing 
and  the  6-cm  disk  diameter.  The  steel  disk  had  a  10^-S/m  conductivity  and  its 
permittivity  was  that  of  free  space,  and  the  remainder  of  Figure  1  had  free-space 
properties. 

Fields  were  computed  with  an  fdtd  program  that  was  written  [6]  in  accor¬ 
dance  with  an  early  manuscript  version  of  [7],  and  whose  absorbing  boimdary 
condition  was  later  replaced  [6]  with  a  Berenger  PML  [8].  This  FDTD  program  is 
second-order  accurate  in  space  and  time.  The  program  was  run  with  CFL  =  1/2 
and  for  several  grid  refinements.  The  discontinuity-related  noise  in  the  Kyi  com¬ 
putation  decreased  as  the  grid  was  refined,  and  the  spectra  of  the  noise  and 
the  main  physical  features  of  the  computed  signals  were  well  separated  on  the 
finest  spatial  grid  (2840  x  5232),  which  yielded  4735  points  per  carrier-frequency 
period.  The  same  2840  x  5232  grid  was  used  to  compute  each  one-cycle-sinusoid 
response  Rj{t)  in  Figure  2  and  to  compute  /iTh- 
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Figure  2.  The  total  signal  received  at  every  other  location  in  Figure  1. 
These  16  signals  have  different  vertical  offsets,  with  the  largest  and 
smallest  offsets  used,  respectively,  for  the  top  and  next-to-bottom  x- 
marked  locations  in  Figure  1.  Note  the  disk’s  shadow  and  the  delayed- 
wave  components. 
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The  integral  kernel  Ks.  and  each  received  signal  Rj  were  both  convolved 
with  the  40-point  filter  as  described  in  Eqs.  (9)  and  (Al).  The  filtered  func¬ 
tions  $  ♦  Rj  and  $  *  were  truncated  as  in  Eqs.  (10)  and  (A7),  and  then 

the  time  origin  was  shifted  so  that  t  =  0  represented  the  pulses’  analytically 
determined  arrival  times  (ta)-  The  filtered,  truncated,  time-shifted  pulses  were 
fed  into  a  second-order-accurate  routine  (Eq.  (5)).  Each  received  signal  was 
deconvolved  separately  from  the  31  other  signals,  resulting  in  32  independent 
inferences  of  the  incident  field  /inc-  A  typical  reconstruction  is  shown  in  Fig¬ 
ure  3,  which  corresponds  to  the  signal  received  at  the  15‘^  cross  from  the  bottom 
of  Figure  1.  Figure  4  shows  that  each  of  the  32  second-order-accurate  recon¬ 
structions  of  /inc  (dash-dotted  curve  with  circles)  reproduces  the  Loo  norm  of 
the  one-cycle  sinusoid  (flat,  boldfaced  curve)  /inc  within  2%.  The  L2  norm  (not 
graphed),  whose  square  is  proportional  to  the  energy  of  a  free-space  pulse,  is 
reproduced  to  within  1.3%  error. 


Figure  3.  Reconstruction  of  the  incident  signal  /nc?  using  the  signal 
received  at  the  IS'^^  location  from  the  bottom  of  Figure  1.  The  incident 
signal  and  its  reconstruction  almost  overlap. 
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Location 


Figure  4.  Restoration  of  the  Loo  norm  for  the  32  locations  in  Fig¬ 
ure  1.  Note  the  shadow  and  diffractive  foci  of  the  received  signals.  The 
reconstructed  signals’  norms  were  within  2%  of  ||/inc||Loo>  which  was  1. 

The  predominant  error  in  these  0(/i^)  results  is  barely  evident  in  Figures  3- 
4.  To  find  the  predominant  reconstruction  error  in  Figure  3,  look  along  the 
left-hand  vertical  axis  for  a  brief,  diagonal  segment  near  0.15  V/m.  Similar 
features  in  all  32  of  the  0(/i^)  reconstructions  predominate  the  relative  errors 
1 1  /reconstructed  -/inc  ||/||/inc||  that  are  illustrated  in  Fi^re  5.  Richardson  extrap¬ 
olation  was  used  [2]  to  obtain  the  O(k^)  results  in  Figure  5.  We  see  there  that 
the  extrapolation  reduced  the  Loo  relative  error  from  ?sll%  to  «8%,  but  that 
the  L2  error  remained  ~  4%.  A  detailed  inspection  of  each  0(h‘^)  time  trace 
showed  that  it  differed  from  the  corresponding  O(h^)  graph  mainly  in  the  first 
extrapolated  data  point.  It  appears  unlikely,  however,  that  these  small,  brief  re¬ 
construction  errors  would  have  any  practical  significance  for  the  radar  problem 
that  these  computations  model.  That  model  radar  problem,  after  all,  has  a  steel 
pipe  (Figure  1)  that  causes  a  deep  shadow  and  diffractive  foci  in  the  received 
signals  (Figure  4),  which  both  the  O(k^)  and  O(h^)  deconvolutions  have  largely 
removed  (Figures  3-5). 
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Figure  5.  Lj,  relative  errors  in  the  reconstructed  signals.  Richardson 
extrapolation  decreased  the  Loo  error  considerably,  but  left  the  L2  error 
essentially  unchanged. 


5.  CONCLUSION 

This  report  is  related  to  [4],  which  used  FDTD  to  compute  the  Kronecker-delta- 
function  response,  which  [4]  then  used  to  propagate  other  incident  fields  _finc 
using  Eq.  (2a).  Reference  [4]  considered  other  reference  pulses  as  alternatives  to 
Ks,  but  it  concluded  that  Ks  was  the  simplest  reference  pulse  to  use.  Thus,  Ks 
and  Eq.  (2a)  are  best  for  the  (convolutional)  forward-propagation  problem  of  [4], 
and  Kn  and  Eq.  (2b)  are  best  for  the  (deconvolutional)  reconstruction  problem 
in  the  present  report.  These  two  combined  methods  Eqs.  (2a)-(2b)  are  useful 
despite  their  integral  kernels  being  FDTD-computed  responses  to  discontinuous 
fields.  The  t- domain-deconvolution  method  of  Eq.  (2b),  for  example,  is  stable 
and  0{h^)  accurate  —  it  is  0{h^)  accurate  after  Richardson  extrapolation  — 
and  it  consequently  may  compete  with  frequency-domain  methods. 
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Significant  work  was  required  because  of  the  numerical  error  in  a  computed 
Heaviside  response  Ks,.  Because  a  discontinuous  integral  kernel  tends  to  con¬ 
serve  significant  digits,  as  shown  by  Eq.  (7),  it  was  best  to  trimcate  (T)  the 
FDTD-computed  K-^  and  then  filter  out  (#)  the  large-ampHtude,  high-frequency 
noise  that  would  otherwise  have  obliterated  the  T-restored  discontinuity.  Those 
steps  wotild  have  been  unnecessary  were  Kn  known  exactly. 

We  finally  consider  what  one  could  do  in  practice,  using  laboratory  data. 
Heaviside-step  pulses  are  problematical  in  the  laboratory,  so  it  is  likely  that  Kb. 
would  have  to  be  inferred  from  other  measurements.  For  removing  interference 
in  an  operating  band  [wmin,  t<^MAx],  for  example,  one  could  presumably  infer 
Kb  by  propagating  a  reference  pulse  /ref  whose  spectrum  covers  [wmin,  ^^MAx] 
and  whose  risetime  is  <C27r/uj.  To  infer  Kb  one  would  then  solve  Eq.  (2b)  for 
Kb  using  the  received  signals  R  resulting  from  the  known  /nc  =  /ref-  If  such  a 
Kb  were  truncated,  then  the  result  could  presumably  be  used  to  reduce  inter¬ 
ference  in  the  operating  band  [wmin?  t<^MAx]-  There  are  two  final  observations: 
The  response  to  a  short-risetime  /ref  can  be  synthesized  by  superimposing  the 
responses  to  several  band-limited  signals.  Second,  as  long  as  the  propagation  re¬ 
mains  linear,  the  techniques  developed  here  apply  regardless  of  whether  multiple 
scattering  is  evident  in  that  linear  propagation. 
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APPENDIX  (filter  and  truncation) 


This  appendix  derives  the  identities  Eqs.  (A4)-(A5),  which  are  used  in  the  devel¬ 
opment  Eq.  (9)  of  this  report’s  numerical  method.  That  numerical  method  also 
uses  a  truncation-related  approximation  that  is  studied  here  in  Eqs.  (A7)-(A16) 
and  shown,  under  a  reasonable  assiimption,  to  vanish  as  the  stepsize  h  tends  to 
0.  The  assumption  is  that  the  superluminal  parts  of  the  FDTD-computed  R  and 
JTh  tend  to  0  pointwise  as  h  — ^  0. 

The  numerical  example  in  Section  4  uses  a  centered  filter 


1 

2m  +  tm 


(Al) 


for  which  differs  from  Tm  by  at  most  one  timestep  (h),  and  min(tni,  2m)  >  0. 
That  filter  is  also  a  convolution  on  (— oo,  oo) 


/OO 

—  s)f{s)ds  (A2) 

-OO 

^(t)  ==  {Tm  +  tm)  ^  H(t -|- 2M)H(im  ~  i)?  (-^3) 

for  a  generic  function  /.  We  also  note  that  /inc(t)  =  /inc(t)H(0)  2?(t) 
R{t)E{t-tc),  and  Knit)  =  KH{t)E{t-Q,  where  U  =  iaCFL  <  is  the  (super- 
Itiminal)  FDTD-computed  pulse-arrival  time  and  ta  is  the  analytical  (lightspeed) 
pulse-arrival  time.  We  will  assume  that  ic  >  2m,  as  is  true  for  the  computation 
in  Section  4.  Next  we  will  show  that 


^*{d;^R)=d-^{^*R)  (A4) 

$  *  (ATh  *  Zinc)  =  ($  *  Kb.)  *  Zinc  (A5) 

identically  for  all  t  G  (— oo,  oo). 

To  derive  Eq.  (A4),  use  Eq.  (Al)  and  the  obvious  definition  of  the  an¬ 
tiderivative  and  use  a  change  of  variables,  to  establish 

($  *  A)  =  T^r—r  /  /  dsJl(s),  (A6) 

I  J —tm  J  s* 
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which  is  identifiable  as  Eq.  (A4)  because  Tm  ^  tc  R{t)  —  R{t)H(t  —  tc)- 
Equation  (A5)  merely  asserts  the  elementary  associative  property  of  convolution. 
The  identities  (Eqs.  (2b)  and  (A4)-(A5))  3deld  Eq.  (9). 

Although  Eqs.  (A4)-(A5)  are  identities  for  any  linear  system,  such  as  phys¬ 
ical  scattering  or  its  FDTD  approximation,  the  computations  in  Section  4  used 
tnmcated  quantities 


To/(t)  =  /(t)H(f-4)  (A7) 

for  /  €  {Kji,  R}  as  in  Eq.  (10).  Such  truncation  introduces  errors 


Ar  =  ar'  [To  ($  *  J?)]  -  *  R)  (AS) 

=  [To  ($  *  ATh)]  *  fine  —  (A9) 

because  the  FDTD-computed  functions  A'h  aiid  R  are  supported  on  [ic,oo]  and 
because  4  =  taCFL  <  <a-  We  will  now  estimate  Ar  and  Aji^,  starting  with  Ar. 

Equation  (Al)  implies 


$*/  =  (#*  /)H(t  4-  Tm  -  tc)  (AlO) 

for  /  €  {R,  Ah}-  Evaluate  Eq.  (AS)  using  Eqs.  (A7)  and  (AlO)  to  show  that 
Ar  =  —  Jtc-TM  *  R)^s,  which  is  independent  of  t.  Apply  Eq.  (Al)  and  a 
change  of  variables  to  the  previous  equation  to  establish 


Tm 


Aii  = 


Tm  + 


l~  tjn  J — J tc  — 


-Tm4-s' 


R(s)dsds'. 


(All) 


Use  Eq.  (All)  and  the  triangle  inequality  to  estimate  lAii|  by  bringing  the 
absolute  value  through  the  integrals,  then  integrating  the  result  over  the  union 
of  all  s'-dependent  intervals  of  integration  in  the  inner  integral  of  Eq.  (All), 
and  then  using  R{t)  =  A(t)H(t  —  tc)  to  obtain 


IAhI  < 


R 


Li  (^cj^a+TM) 


(A12) 


The  previous  norm  is  over  an  interval  that  is  only  slightly  larger,  by  an  amount 
Tm,  than  the  support  of  the  superluminal  part  of  the  FDTD-computed  R.  The 
assumption  in  this  appendix’s  first  paragraph  and  the  observation  that  Tm  — >  0 
as  h  — ^  0  (because  Tm  =  20h  in  Section  4)  then  imply  that  Ar  0  as  h  0. 
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Now  we  estimate  Ak  in  the  h  — s-  0  limit.  Use  the  Young  theorem  [9]  and 
the  identity  To  ($  *  Ks)  =  H(t +  Tm  — <c)H(t — <a)($  *  K-r)  to  evaluate  Eq.  (A9) 
and  obtain 


< 

#*a:h 

fine 

Ls(-C30,  00) 

haltc^TMi  ^a] 

Lp[0,  <a] 


(A13) 


for  all  s,  p  €  [1,  oo],  where  a  ^=1+3  ^  —p  We  will  use  Eq.  (A13)  to  show 
limft-^o  i^K  =  0  by  establishing  that  there  is  a  cr  for  which 


lim 

h-*Q 


$  *  ATh 


=  0. 

L(t[<c“7mi  *a] 


(A14) 


Indeed,  we  will  show  that  Eq.  (A14)  applies  for  all  cr  G  [1,  oo),  the  exclusion  of 
C50  being  notable.  For  <7  =  00  one  can  proceed  in  the  most  straightforward  way, 
using  Eq.  (Al)  and  =  JiCh(0H(^  ~  to  obtain 


$  *  ATjfi 


< 

Loo[“,  i] 


[a— tm,  6+Tm] 


(AI6) 


for  any  a  <  b.  The  cr  =  00  case  is  excluded  from  Eq.  (A14)  because  Eq.  (A15) 
would  bound  the  Loo[tc  —  Tm,  ta]  norm  of  #  *  ATh  from  above  by  the  Loo[tcj  + 
Tm]  norm  of  Ks,  which  cannot  be  forced  to  0  as  h  — >  0  because  of  the  h- 
independent  discontinuity  at  i  =  ta  of  the  exact,  analytical  Kji',  but  Eq.  (A15) 
will  soon  be  used  as  a  lemma.  Because  Ke  is  discontinuous  at  <  =  <a,  we  now 
evaluate  the  L<,^(tc  —  Tm,  ta)  norm  (cr  <  cxj)  of  $  *  ATh  by  estimating  separately 
the  ^u-Tm  components  of  the  definition  of  that  Ler  norm.  The 

result,  after  some  straightforward  manipulation  involving  Eq.  (A15),  is 


$  *  ATh 


(7 

Lcrftc—^M,  ^a] 


<(ta  -  tc) 


Kn 


+ 

Loo[<c ,  ^a] 


Tm 


cr 

Loo[ia~?M“-tm,  ^a+TM] 


(AI6) 


for  all  a  £  [1,  co).  The  first  term  in  the  sum  in  Eq.  (A16)  vanishes  as  h  — »•  0 
because  of  the  assumption  in  this  appendix’s  first  paragraph.  Then  lim/i_^o  Tm  = 
0  implies  that  the  left-hand  side  of  Eq.  (A16)  vanishes  in  the  same  limit;  and 
Ak  0,  too,  because  of  Eq.  (A13).  We  have  now  shown,  subject  to  the 
assumption  at  the  beginning  of  this  paragraph,  that  max(|AK|,  |Ajtc|)  — >  0 
as  h  0.  Thus,  the  error  introduced  by  truncation  in  Eqs.  (10)  and  (A7)- 
(A9)  vanishes  as  h  — >  0,  provided  the  underlying  fdtd  method  is  reasonable,  as 
defined  in  the  first  paragraph  of  this  Appendix. 
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